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Abstract 

The elastic constants of solid hep 4 He are investigated in molar volumes ranging from about 
melting up to approximately 14 MPa. Properties of interest are determined by averages formed 



T3 

a 

q ■ from computed values in configurations sampled by Monte Carlo of a model wave function. Devi- 



ations from known elastic relations are reported near the density where the supersolid fraction is 
maximum. The results offer further evidence that the supersolid state is related to elastic constants 
anomalies and that both are manifestations of a single physical process. 



For more than a hundred years, since it has been liquefied, helium has shown intriguing 
behaviors and has provided invaluable clues for the condensed matter physics. Nowadays, 
despite a reasonable good understanding of its liquid phase, solid helium is still revealing 
fascinating aspects that defy experimentalists and theoreticians. The shear modulus soft- 
ening with temperature increase, as observed by Day and Beamishfl], and its connections 
to supersolidity, superfluidity superimposed to crystalline order, as seen in the pioneering 
experiments by Kim and Chan [2, 3] illustrate very well this situation. 

Elasticity is ubiquitous among materials and has been extensively investigated in solids 
where atomic exchange is negligible, and the atoms are reasonable localized and distinguish- 
able. If forces that are not too strong are applied in the surfaces of these solids, they will 
produce mechanical deformations that will increase the potential energy. As the solid is 
released from these forces and assuming they are small enough, it will go back to its original 
state. However, what can be said about quantum solids? Systems where the zero-point 
kinetic energy is comparable to the binding energy of the atoms that perform large excur- 
sions from their mean positions in the lattice. How will these quantum systems respond to 
applied forces? 

This letter reports values of all the five independent elastic constants, as a function of the 
density, which characterize a defect-free single crystal of 4 He in an hep structure. We have 
analyzed our results using relations among the elastic constants that are satisfied under well 
defined conditions. Since all the elastic constants were determined, it was possible to some 
extent verify the internal consistence of our calculations and search for elastic anomalies 
in this quantum solid. In fact, we have observed deviations of the relations occurring near 
the same density where the mass decoupling in the movement of torsional oscillators has 
its maximum. [4] In our view these results are further evidences that supesolidity[2, 3] and 
elastic anomalies [1, 5] of the shear modulus are related phenomena and consequence of a 
single physical process. 

The elastic constants were calculated at zero temperature using a model wave function 
that has allowed an excellent description of various properties of 4 He in both the liquid and 
solid phases. [6, 7] It is worth mentioning that the interacting potential for the helium atoms 
is very accurately known as more than one study has shown. [8, 9] These facts, together with 
Monte Carlo calculations performed with a reasonable number of atoms, assure that the 
estimated quantities can have good reliability. 



A system formed from 4 He atoms of mass m can be well described by the Hamiltonian 
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where V(r) is the Hartree-Fock dispersion two-body (HFD-HE2) potential of Aziz and col- 
laborators. [10] This inter-atomic interaction is still widely used, [11, 12] even if more accurate 
potentials exist. [13, 14] The HFD-HE2 potential was kept in our calculations because it does 
not seem that our results will be significantly modified if the more recent potentials were con- 
sidered. Moreover there is already a large body of work performed with the model function 
used in this work that relays on the HFD-HE2 inter-atomic potential. 

We characterize the solid phase by a shadow wave function[15, 16] (SWF) written as 



N 



V{R) = ]Jexp 



%<j 



1 (b^ 



I \ Tij 



r N N 

/ dS JJexp [-C\ri - s,| 2 ] JJexp [-8V(a Sij )\ (2) 



Kj 

where R = {r l5 r 2 , . . . , r^}, S = {si, s 2 , . . . , Sjv}, and dS = {ds ly ds 2 , ■ ■ ■ , ds N }. The in- 
tegrations over the auxiliary variables, also called shadow particles, are performed in the 
whole space. The parameters b, C, 5 and a are chosen so that the expectation value of the 
energy is minimized. The correlations between the shadow particles are imposed through 
the HFD-HE2 potential rescaled in its amplitude and inter-particle distance through the 
parameters 5 and a, respectively. This model wave function is considered one of the best 
available descriptions for systems consisting of helium atoms. [6, 7] It has been used to study 
several helium properties[12, 17, 18] including the structure and mobility of linear defects in 
a calculation of the Peierls-Nabarro stress. 

Elastic constants are computed using the Parrinello and Rahman method. [19, 20] This 
treatment allows periodic boundary conditions to be considered together with distortions 
that must be applied in the simulation cell in a easy way. If the orthogonal vectors L x , ~L y 
and L 2 define this box along the directions (x,y,z), a position r, in the cell can be written 
as 

I-; = &L X + r]iL y + dL z = h Si , (3) 

where the coefficients £j, r\i and Q range between and 1 and can be considered components 
of a vector Sj. Then the vector rj can be written as a product of a matrix h by s$. The 
matrix li is diagonal with elements equal to the modulus of the vectors that define the box. 



Homogeneous deformations change the positions iv In our calculations the new positions 
are obtained through a non-diagonal matrix li 

r[ = lis, = Mio" 1 !-, . (4) 

The strain-tensor then can be written as 

*<# = \ (14- W - 1) , (5) 

where the Greek symbols stand for any of the coordinate axes, Q = Wh. denotes a metric 
tensor and I is the unit matrix. 

The stress-tensor a a p is calculated through the virial theorem. [21] It is estimated by 
averages of values calculated from configurations sampled from the square of the SWF using 
the Metropolis algorithm. In the linear elastic regime, the stress and strain tensors are 
related by the generalized Hooke's law 
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where the c a p lx are the stress-strain coefficients. The elastic constants, 
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were obtained by computing the slopes of stresses as a function of strains for selected defor- 
mations that have a single strain-stress coefficient different from zero. Elastic constants are 
equal to the stress-strain coefficients if no hydrostatic pressure P is applied to the crystal. 
Otherwise the elastic constants are related to stress-strain coefficients by[22] 

Relations among the elastic constants were used [23] in the past to complete their individ- 
ual determination from experimental data. Here, since we are able to calculate all of these 
constants independently one from the others and with good accuracy, the relations are used 
to validate the internal consistence of our calculations and to look for deviations that might 
signalize elastic anomalies. 

The generalization of Hooke's law considered together with the symmetric character of 
the stress tensor, the reversibility of the deformation work and the symmetry of a hep crystal 



result in a system with only five independent elastic constants. If the additional constraint 
of a short-range central interacting potential is assumed, the number of elastic constants is 
reduced to four independent values. This additional symmetry is one of the Cauchy relations 
that is still valid in an hep structure [24] 

C13 = c 44 . (9) 

From now on we use the Voigt notation. [25] It is worthwhile to note, as it has been extensively 
discussed in the literature, that the Cauchy relations are not rigorously fulfilled by any 
material. Nevertheless this fact has not prevented its use as a source of informations, for 
instance, about bonding properties of the system. [23, 26] 

Another relation among elastic constants can be obtained by considering the linear com- 
pressibility. It is given by the relative longitudinal change in length of a given direction 
with respect to the hexagonal axis when the system is under unit hydrostatic pressure. In 
general this quantity is non-isotropic. However if the ratio c/a of the lattice parameters 
along the hexagonal axis c and the basal plane a is independent of the pressure, the linear 
compressibility is isotropic and as a consequence the relation 

en + C12 = c 33 + c 13 (10) 

is satisfied. 

We have computed properties of interest at several molar volumes resulting in densities 
varying between 0.0292 and 0.0353 A -3 . Simulation cells we have used contained 720 atoms 
in an ideal hep structure. The equilibration process consisted of 2 x 10 4 Monte Carlo 
(MC) sweeps (attempts to move randomly all the particles and shadows coordinates). At 
all densities, the parameters of the SWF are those that minimize the expectation value of 
the Hamiltonian operator of Eq. (1). The variational search in the parameter space was 
performed either by us or is reported in the literature. [16] 

The stress tensor components were estimated by performing simulations of at least 1.6 x 
10 6 MC sweeps for each of the different strains. Each elastic constant was determined 
using four different distortions in the simulation cell. The maximum strain levels we have 
considered was 0.25%. 

The elastic constants were determined in molar volumes ranging from 20.66 cm 3 , near 



the melting density, to 17.06 cm 3 . The results are displayed in Fig. (1) together with the 
data from the literature. [23, 27-31] The general agreement between theory and the few 
experimental data points is good. From the figure we can see that constants Cn, Cyi and 
C13 have variations with respect to density that clearly show some sort of discontinuity near 
the molar volume 19.36 cm 3 . Although constants C33 and C44 seem also to show this feature 
at the same density, it is much less pronounced. It is worthwhile mentioning that near this 
density, a maximum in the supersolid fraction [4] has been observed. 

After completing this work we became aware of Ref. 6 where elastic constants were 
investigated by diffusion Monte Carlo. In this paper the authors point out that the elastic 
constants that we have already published[12] at a single density p = 0.0294 A -3 is in 
good agreement with their results. At other densities we see a less satisfactory agreement. 
However the most important disagreement between our results are the overall behavior of 
the elastic constants with pressure. While in Ref. 6 the elastic constants increase linearly 
with pressure, our results are far from this behavior. Especially for the elastic constants Cn, 
C\i and C13, we see roughly constants values up to p = 0.0311 A~ 3 where the magnitudes 
start to increase with a character that is hardly linear. We do not know the reasons for such 
discrepancies. However it would be difficult to attribute them to our model. It has been 
recognized as excellent. [6] Eventually they might be due to finite-size effect. Although the 
authors of Ref. 6 have made a careful analises of this effect on the energy, it is well known 
that simulations made to obtain the elastic behavior of materials require large systems. Our 
system is more than three times bigger than theirs. 

In Table I at each one of the given densities, we report values determined in our calcu- 
lations for the elastic constants, the total energies per atom and the pressures, obtained by 
taking the trace of the stress-tensor. 

We have analyzed our results considering the relations among the elastic constants already 
discussed. For this purpose it is useful to rewrite Eq. (9) using Eq. (8) and define 

(C u + P)- (C 13 - P) . . 

s = few) ' (11) 

this is a quantity that depends on the shear modulus jjl = C44. In principle 5 is zero if 
zero-point motion effects and many-body components in the inter-atomic potential were 
negligible. As already point out, the general validity of the Cauchy relation that originates 
Eq. (11) might be questionable. However whatever value 5 could have it would be hard to 
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1. Elastic constants of 4 He (solid symbols) as a function of the molar volume in a log-log 
The empty symbols stand for the experimental data [23, 27-31]. Values are from the left or 
hand scales according to the arrows in the figures. The lines are guide to the eyes. 



expect abrupt changes in its magnitude as the pressure is slightly varied. Fig. (2) presents 
the values we have computed for 5. First of all, let us note that at both ends of the densities 
range we have considered, 5 is approximately 0.4. A possible value, because in a system 
made from 4 He atoms the requirements for the validity of the Cauchy relation, 5 = 0, would 
hardly be met. However near the molar volume 19.36 cm 3 , 5 has a peak with a maximum 
value equal to almost 2.0. This is near the density where the supersolid fraction has a 
maximum value. [4] The quantity § has also been computed using the available experimental 
data. Values obtained from experiment are in good agreement with theory. However there 



p (A -3 ) V (cm 3 ) E (K) P (MPa) C u C 33 C* 44 C 12 C 13 

0.0292 20.62 -5.0770 ± 0.0010 2.4529 ± 0.0047 59.3 ± 1.2 78.1 ± 1.9 16.63 ± 0.56 29.6 ±1.5 14.89 ± 0.68 

0.0294 b 20.48 -5.0349 ± 0.0007 3.0957 ± 0.0026 60.8 ± 1.0 77.9 ± 3.4 17.10 ± 0.34 34.4 ± 2.2 14.40 ± 0.85 

0.0301 20.01 -4.8646 ± 0.0007 4.0707 ± 0.0028 65.7 ± 1.8 82.9 ± 2.3 18.07 ± 0.34 32.9 ±3.1 15.8 ± 1.3 

0.0307 19.62 -4.6410 ± 0.0007 5.0941 ± 0.0029 65.5 ± 1.2 93.5 ± 2.2 18.83 ± 0.22 32.1 ± 1.9 14.81 ± 0.67 

0.0311 19.36 -4.4799 ± 0.0008 5.7905 ± 0.0031 63.6 ±1.8 93.1 ±5.5 19.43 ± 0.39 31.7 ±2.0 14.7 ±1.1 

0.0316 19.06 -4.2864 ± 0.0007 6.5933 ± 0.0030 82.7 ±1.0 105.9 ± 3.7 22.38 ± 0.36 41.9 ± 1.6 21.3 ± 1.9 

0.0320 18.82 -4.0713 ± 0.0007 7.2055 ± 0.0031 95.4 ± 2.6 116.9 ± 1.0 24.29 ± 0.39 50.9 ± 1.2 28.61 ± 0.68 

0.0329 18.30 -3.5904 ±0.0017 8.8184 ±0.0081 117.0 ±4.2 142.1 ± 3.7 27.10 ± 0.35 65.1 ± 3.4 36.3 ± 1.5 

0.0341 17.66 -2.8139 ±0.0006 11.6933 ±0.0035 148.0 ±2.5 176.3 ± 1.3 35.77 ± 0.22 78.9 ± 1.2 48.16 ±0.95 

0.0353 17.06 -1.9265 ±0.0013 14.1040 ±0.0074 171.3 ±4.3 202.4 ± 8.7 41.82 ± 0.57 90.4 ± 2.7 52.8 ± 2.7 

b Results from Ref. 12 

TABLE I. Elastic constants Cij in MPa calculated at the given densities (first column) in solid hep 
4 He. The second column shows molar volumes. The third and fourth columns show the values of 
energies per atom and pressures. 

are experimental data only at few densities, what precludes a detailed comparison between 
theory and experiment at the density where 5 show its peak. 

The behavior of 5 can not be imputed to an artifact of our model. As a possible mean 
to verify this fact we have considered the equation of state (EOS) of the total energy per 
atom as a function of the density. The computed EOS is in very good agreement with 
experimental data as it is possible to verify in the inset of Fig. (2). The computed pressures 
are within the values one expects for this quantity. It is also displayed in the inset, and it 
does not show neither any signal pointing towards an specific difficulty of the model. 

For solid 4 He as we can learn from experiment [33] and theory, [6] the axis ratio c/a is 
almost pressure independent and very near 1.63, the value c/a has in an ideal hep crystal. 
In this context, we expect the relation expressed in Eq. (10) to hold. Because of Eq. (8), 
this relation is valid for both the stress-strain coefficients and the elastic constants. We have 
computed the ratio 1Z = (Cu + C , i 2 )/(C , 3 3 + Ci 3 ), as a way of verifying the validity of the 
relation of Eq. (10) and plotted our results in Fig. (3). Since 4 He has the ratio c/a practically 
constant and independently of the pressure, we expect the ratio TZ to have also a constant 
value equal to one. In fact we see this behavior at most of the densities where the calculations 
were performed. However, as we see in the figure, TZ deviates from this value. This quantity 
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FIG. 2. The ratio 5 that assess deviations of the Cauchy relation as a function of the molar volume 
and density, respectively lower and upper scales. The empty circles stand for the experimental 
data [31]. The inset shows the computed energies per atom as a function of the density (solid 
triangles and left-hand side scale) and the experimental data [32] (open symbols). The line is 
a fourth order polynomial fit to the energies. The pressure (solid squares) has its units in the 
right-hand side. The line is obtained through the derivative of the EOS, P(p) = p 2 d(E)/dp 



has a minimum in a form that resembles a peak upside down at the same density where the 
supesolid fraction and 5 of Eq. (11) have their maxima. We also have computed the ratio 
TZ using experimental data and once again only in a few densities there are available data. 
This is not a helpful situation for a direct comparison between experiment and theory where 
we see the largest deviation of 1Z. However where the comparison is possible the agreement 
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FIG. 3. The elastic constants ratio (Cn + C\2)/{C^ + C13) as a function of the molar volume 
at the lower scale and the density in the upper scale. The solid symbols stand for our theoretical 
results and those empty for the experimental data. [31] 

is good. 

In this work, as we have described a hep defect-free single crystal of 4 He, the elastic 
constants were computed as a function of the density. The deviations from the relation 
of Eq. (10), derived from the well established fact that 4 He has an isotropic linear com- 
pressibility, show according our view that our calculations are consistent with experimental 
observations that a system formed from these atoms present elastic anomalies. At the same 
time this relation being obeyed at the lowest and highest densities we have considered offer 
a confirmation of the consistence of our calculations. The behavior 5 defined in Eq. (11) 
corroborates this fact. The values of this quantity are able to display the zero-point mo- 
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tion effects at the extremes of the density range and an elastic anomaly, as exhibited by 
a relatively large peak. We do not believe that values of 5 ~ 0.4 could be attributed to 
a non-central character of the atomic interaction. Atomic three-body interactions can be 
accurately computed [34] and we know [9] they are small. It would be desirable to have addi- 
tional experimental data for the elastic constants near the density p = 19.36 A~ 3 to confirm 
our finds. 

Our results show that anomalies of the elastic constants are properties of a defect-free 
single crystal of 4 He in an hep structure. We stress that the calculations have not considered 
3 He atoms impurities or the presence of defects. The relations among the elastic constants 
we have used present deviations consistent with the behavior of the supersolid fraction as a 
function of the system density, an increase of its value up to near p = 0.0311 A -3 followed 
by a decrease in its value. These facts corroborate that supersolidity and shear anomalies 
are related phenomena probably a manifestation of a single physical process. 

Comprehension of elastic properties of solids are important for materials engineering and 
for basic science as well. Elasticity is one of the next chapters where the investigation of 
solid 4 He will bring important contributions to quantum many-body and condensed matter 
physics. 
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high-performance computing facility at Universidade Estadual de Campinas and at the lab- 
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